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Abstract 

Starting from the Pauli Hamiltonian operator, we derive a scalar quantum kinetic equations for spin- 1/2 systems. 
Here the regular Wigner two-state matrix is replaced by a scalar distribution function in extended phase space. Apart 
from being a formulation of principal interest, such scalar quantum kinetic equation makes the comparison to classi- 
cal kinetic theory straightforward, and lends itself naturally to currently available numerical Vlasov and Boltzmann 
schemes. Moreover, while the quasi-distribution is a Wigner function in regular phase space, it is given by a Q- 
function in spin space. As such, nonlinear and dynamical quantum plasma problems are readily handled. Moreover, 
the issue of gauge invariance is treated. Applications (e.g. ultra-dense laser compressed targets and their diagnostics), 
possible extensions, and future improvements of the presented quantum statistical model are discussed. 

PACS: 52.25.Dg, 51.60.+a, 71.10.Ca 
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1. Introduction 



Quantum kinetic theory has a long history. In many respects, it all started with the seminal paper by Wigner in 1 932 
|Q]] , see also Refs. 12, 0] , and the later developments of Moyal jH] • While the approach of Wigner has the advantage of 
being of interest for the interpretation of quantum mechanics Q3Q , and also for the development of quantum optics (for 
an overview, see e.g. JfJ), detailed calculations of material properties in condensed matter systems have relied to a 
large extent on either semiclassical techniques [7], in which the collisional operator in Boltzmann's equations involves 



as well as diagrammatic techniques CI. The 
121 1. has been successful in dealing with certain 



quantum transition probabilities, or Green's function techniques 
theory of Baym and Kadanoff, as well as the works of Keldysh yj , 
quantum transport phenomena. The theory contains memory effects (nonlocal terms, both in space and time), has a 
straightforward interpretation in terms of the different Green's functions, and the theory works well even on time- 
scales shorter than the typical relaxation time of the system in question. However, the gap between classical plasma 
physics and quantum transport theory does not seem to have been bridged, probably due to reasons of formalism as 
well as a difference in application of the respective models. Moreover, while the Kadanoff-Baym equation gives a 
very good description of certain systems, it is perhaps not well-suited to some of the future applications of quantum 
kinetic theories, such as high intensity la ser-p lasma interactions [13], high energy density physics [ 14], and nonlinear 
collective quantum problems 1 15 , Kjl 

In particular, the field of quantum plasmas has recently attracted, a perhaps unexpected, interest in the field of 
laser plasmas [20, 21, 22], where high density ionized plasmas can be created in the laboratory. Moreover, the event 
of nano-devices and technology on sub-micron scales, such as quantum dots f23, 24, 25t] and plasmonic components 
1 261, 27ll . has sparked the interest of many researchers of analyzing the dynamic and nonlinear properties of such 
systems. A recent result is that quantum effects in plasmas can be important in parameter regimes that for a long time 
have been considered purely classical B28I1 . 
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The above discussion is mainly related to the statistical and dispersive behavior of unmagnetized quantum plasmas 
02911 . However, one intrinsic non-classical property of quantum systems is the spin. The magnetization that follows 
from the intrinsic spin, as well as that of orbital angular momentum, is of course the foundation for many important 
material properties J30ll . Investigations of such condensed matter systems are often directed towards equilibrium 
properties, although the nonlinear dynamics of magnetization is sometimes interest and probed using the Landau- 
Lifshitz-Gilbert equation [31]. There are a variety of different physical systems where the spin can be of importance, 
such as metal alloys and semiconductors material for memory use [31], cold atom gases IB2I1 . and high density and 
high field astrophysical plasmas ifTill . to mention a few. Collective effects originating in the plasma particle species 
spin has therefore recently become an active field of research for fully ionized systems (see e.g. B19LI33I1 and references 
therein), in particular in the nonlinear regime, where spin solitons IB411 and ferromagnetic behavior in plasmas can be 
found [35]. Many of the studies presented in the literature have so far been of a theoretical nature, but it is not difficult 
to envision future applications to e.g. plasmonic devices 12611 or femtosecond physics Pol . 

For the purpose of connecting classical plasma physics to the evolution of nonequilibrium quantum systems, 
utilization of quasi-distributions is of great value. First of all, the interpretation of the quasi-distribution function 
using ensemble averages of observables is in direct analogy with the classical case. It is even possible to directly 
construct a quasi-distribution, such as the Wigner function from measurements f37ll (with the only information loss 
being the initial phase). Second, the quasi-distribution evolution follows from the quantum Liouville equation for 
the density operator, and gives a quantum analog of the Vlasov or Boltzmann equation. This may also render a 
quantum kinetic theory for the quasi-distribution function useful for adaption of classical numerical codes to the 
quantum regime. There are of course infinitely many ways to construct a quasi-distribution function, giving certain 
elementary requirements (see next section). However, a few quasi-distribution functions are more prominent in the 
literature than others. The best known quasi-distribution function is probably the Wigner distribution Jit], but there are 
many others frequently used. In short, different definitions correspond to different operator ordering, so depending on 
the application different definitions are natural. For example, when considering optical coherence normally ordered 
operators occurs naturally and hence the Glauber-Sudarshan P-distribution 11381 13911 is a convenient choice. On the 
other hand, anti-normal ordered operators i.e. the Q-function or the more general Husimi functio n 1I40I1 are useful 
when dealing with quantum chaotic systems. For reviews of the subject see for example Refs. I5 U4U 14211 . 

In this paper, we will construct a quasi-distribution function for a particle with spin- 1/2 as a combination of a 
Wigner distribution for the position and momenta and the Q-function for the spin degree of freedom. Moreover, 
a quantum kinetic equation giving the evolution of this scalar distribution function, in the mean field or Hartree 
approximation, will be derived and applications to magnetized systems will be presented. A discussion of possible 
future applications and research directions will also be given. 

The structure of the paper is as follows. In Sec. 2 we give a short overview of different quantum quasi-distributions. 
In Sec.[3]we consider the evolution equation for a density matrix for a spin-1/2 particle in an external electromagnetic 
field. In Section|4]we go on to derive a combined transformation for the phase space and spin variable. This transfor- 
mation then renders an evolution equation for the system in extended phase space (x, p, s) which is derived in Section 
[5] The extension to the mean field approximation is reviewed in Section [6] and in the following section we calculate 
the thermodynamic equilibrium density matrix for a set of N noninteracting particles. In Section [8] we consider the 
evolution equation in the long scale length limit and compare our results to previous semi-classical kinetic descrip- 
tions in the literature. In Section[9]we consider the linear solutions to the derived equations. Section 10 is devoted to 
a discussion of gauge properties and the fully gauge invariant evolution equation is presented. Finally we summarize 
the main results and discuss future development and applications in SectionfTTI 



2. General requirements of quasi-probability distribution function 

2.1. Historical note 

Following the success of the classical theory of non-equilibrium statistical mechanics, it was natural to seek 
a similar theory for quantum systems in the late 20s and early 30s. However, while the classical Liouville equation 
generates trajectories in phase space as in a classical Hamilton-Jacobi theory, we in the quantum realm have to consider 
the Heisenberg uncertainty principle. This will not allow us to describe, as in classical systems, precise trajectories, 
but rather "smeared out" paths in what would be the corresponding phase space. Indeed, the attempts by de Broglie, 
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Bohm, and others to give a close-to classical interpretation of the Schrodinger equation by Hamilton-Jacobi theory 
shows that, if one is inclined to stick to this interpretational scheme and extend this to statistical interpretations, one 
has to consider the wave function rather as an ensemble of (nonclassical [43, 44]) trajectories (a similar conclusion 
can be drawn from path integral J45ll as well as Ehrenfest techniques ll46l0 . satisfying certain initial and boundary 
conditions. Thus, the introduction by Wigner of a quasi-distribution function (see below) was a natural step in the 
direction of relating measurements to classical transport theory. This is perhaps most obvious in the field of quantum 
optics, where phase space techniques since long has been widely used. Three main definitions of quasi-distributions 
can be found in this field, namely the Wigner function [1], the Husimi (or, equivalently, the Q-) function lf4Qh . and the 
Glauber-Sudarshan P-distribution B38ll39ll . Below we will give a short summary of some of the properties of the first 
two types of quasi-distribution functions (the P-distribution will not be used in the present work). 



2.2. Basic requirements 

Some basic requirements can be imposed on a quantum probability distribution function in phase space, in order 



for it to have a reasonable interpretation B44 14711 . We denote the quantum (quasi-)distribution function by /(x, p) 
(for the moment, we drop the explicit time-dependence for notational convenience) for a given quantum state p of the 
system. Then the marginal distribution functions (x|p|x) and (p|p|p) should be related to /(x, p) according to 

/(x) = j d 3 pf(x,p) = (x\p\x), (1) 

and 

/(p) = j <i 3 x/(x,p) = <p|p|p>, (2) 

respectively. Moreover, we should require that the distribution function is positive definite, i.e. 

/(x,p)>0. (3) 

However, it can be shown that the conditions dTJ — ([3j is not sufficient to uniquely determine a suitable quantum distri- 
bution function in phase space. In fact, Cohen B48I1 has shown that there are infinitely many function /(x, p) satisfying 

tD-0. 

A more complete list of properties that are desirable is found in I42I1 . where expect for the three properties above, 
the additional properties that the distribution function is real, bilinear in the wave function and that the distribution 
functions for eigenstates of the Hamiltonian form a complete and orthogonal set. In fact, it can be shown that in general 
one cannot find a distribution function that satisfies all of (H}-© simultaneously, if one requires the distribution 



function to be bilinear in the wave function [75]. 

Though the above conditions are important when it comes to interpreting the distribution functions a perhaps more 
important condition is that it should be possible to calculate the expectation value of any operator. This condition is 
important since it means that all physically relevant information is included. To calculate the expectation value one first 
map the operator to the corresponding phase space function O = 0(x, p) — > 0(x, p), using the Weyl-correspondence, 
and then calculate the phase space average weighted by the distribution function 

0> = j d 3 xd 3 pf(x,p)0(x,p). (4) 



The mapping from the operator space to phase space depends on which distribution function is used (see Ref. [42 
for a details). Below we will collect the properties of two distribution functions of interest in our context, the Wigner 
distribution [1] and the Husimi function 114011 (or £>-function §Q |47|1 ). These are also perhaps the most frequently 
encountered quantum probability distribution function in the literature (see Ref. 16[] and B47I1 for further references 
and other prominent distributions used in quantum optics, such as the P-distribution of Glauber [38!] and Sudarshan 
j39ll . and their interrelations). 
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2.3. The Wigner function 

The Wigner function for a quantum state p is defined as the Fourier transform of the two-point correlation function 
(i.e., density matrix). Thus, we accordingly have 

P) = TtW f Ae* y/S <x + y/2|p|x - y/2). (5) 
(2nhy J 

Through this definition of the Wigner function, we see that it satisfies the marginal distribution requirements (03 and 
(|2). However, it does not satisfy the positivity criteria (01. The latter property then prevents a probability distribution 
interpretation. However, the negativity of the Wigner function is limited in the sense that the proper number density 
in physical space is n(x) = J d 3 p fw(x, p) which is thus always positive. For a pure state p = \t//)(iff\, this definition 
gives 

fw( x > P) = 7T~~i\3 f dV ip - y/ V (x + y/2)<A(x - y/2). (6) 
(27m) J 

One of the important properties of the Wigner function is that it cannot have too sharp peaks, expressed by 



d\d } p[f W (x,p)] 2 <——r, (7) 

(2nhy 



a result of the noncommutativity between coordinate and momentum operators. 

The time evolution for the Wigner function in an external (analytic) potential V(x, f) is given by 

d ^ + l. V j w W s JH x .j\ fw = , (8) 
at m n \2 J 

where the sin-function is defined in terms of its Taylor expansion in the case of analytic potentials, and we have used 
the indices x and p on the V to denote its operation in phase space. To find the phase space function that corresponds 
to a given operator we must first express the operator in Weyl order ]2[], i.e. express in symmetric products of x, and 
pi, i — 1, 2, 3, using the commutation relations and then substitute jc, — > x and />, — > p. For example, calculating the 
average of the operator Xipj, we have 

1 / \ ih ih 

x iPj = 2 [XiPi + PJ Xi ) + ~J 6i J ~* X P + J 6i J> ( 9 ) 

where 5,j denotes the Rronecker delta function, and hence 

(xiPj) = j d 3 xd 3 p f w {x, p) Lcipj + ^(5 y | . (10) 

2.4. The Husimi function 

The Husimi function (see (fTTT i below) is based on minimum uncertainty wave packets, and it does not satisfy 
(1) and (2) but is positive definite (thus satisfying (3)). As will be seen below, this allows probability distribution 
interpretation of the Husimi function; however, it gives a different greater uncertainty measure than expected through 
naive application of the Heisenberg uncertainty relation. These properties can be immediately understood from the 
following definition. For a given Wigner function, the Husimi function can be obtained through a Gaussian smoothing 
as 

/ H (x, p) = -L f f d\'d?p' exp[-(x' - xf/2d 2 ] exp[-/z 2 (p' - p) 2 /2d 2 ]f w (x' , p'), (11) 
(m) J J 

where the parameter d sets the scale of the smoothing. 

While the Husimi function is positive definite, and produce the correct expectation values of observables, it satis- 
fies an indeterminacy relation of the form 

(Ax) H (Ap) H > h, (12) 
4 



as compared to the relation 

AxAp > H/2 (13) 

for a quantum state (the latter being satisfied by the Wigner function). This results is due to the smoothing introduced 
in the definition of the Husimi function. The Husimi function does not give the probability for the particle to be 
at a certain phase space position, but rather the probability to find the particle in the minimum uncertainty state 
centered around the phase space point in question [49]. Introducing minimum uncertainty states |xo, Po) which satisfies 
Ax 2 Ap 2 = 7z 2 /4 one can write the Husimi function as 

/ H (x,p) = <x,p|p|x,p>. (14) 

However, as mentioned above, it can still be used to calculate any observable, but the operator ordering rule is more 
complicated than in the Wigner case and we will not consider this further here. The evolution of the Husimi equation 
can be found from ([8]) and the definition dTTT) . It is fairly complicated (see 14211 1 and it is more convenient to compute 
the evolution of the Husimi distribution function by evaluating the Wigner function for all times through ([8]). This 
said, we note that although the evolution equation for the Husimi function is more complicated than the corresponding 
equation for the Wigner equation, it is sometimes the convenient choice. One such example is when considering 
chaotic system in which the phase space distribution function becomes ve ry comp licated. The Husimi function, being 
a Gaussian average, may then behave more regularly (see, e.g., Refs. f50l 51. 521). 




2.5. Quasi-distribution functions for spin 

Similarly to the case for phase space it is possible to construct quasi-distribution functions for the spin degree of 
freedom. This has been done already in the 1950's by Stratonovic 15311 . Later on the spin quasi-distribution functions 
were further developed and were applied to problems related to calculating correlation between spins f54l l55l 
The spin quasi-distribution function has also been discussed in connection with quantum scattering problems \51 

As in the case of the regular phase space variables x and p, there is no unique way to introduce a spin quasi- 
distribution function. Scully and Wodkiewicz 15811 give a very good review of the many different choices that can be 
considered. There are at least three different methods for defining spin distribution functions: delta distributions, dis- 
tributions based on coherent states (Q and P) and Stratonovic distribution functions. However, the different outcomes 
of these choices overlap. 

In this paper we will consider only the Q-function for spin which is defined as 

/(0,^) = <s|p|s>, (15) 

where Is) is the state which has spin up in the direction of the unit vector s-s(0, tp) often called a spin-coherent state 
1 591 6011 . Note that this is analogous to the definition of the Q-function in position/momentum space, the latter given 



by Eq. ( TT4"1 >. As for the Husimi or Q-function in phase space, this distribution does not give the correct marginal 
distributions. This means that integrating over the <p angle does not leave the correct distribution function for the 9 
variable. However, it still contains all the information about the system and it can be used to calculate the expectation 
value of any observable, just like in the density matrix formalism. The mapping between spin operators and the 
corresponding spin-space functions will be considered in detail in Section [4] The main reason for choosing to work 
with this particular distribution function for the spin variable is that it is a function on the unit sphere and hence 
resemble the classical picture of a dipole moment, in fact, the evolution equation in the long scale length limit (see 



Eq. d63l ) is almost identical to an equation derived previously from a semiclassical treatment of the spin [61]. The 
Q-function for the spin is also nonnegative which may be desirable in some cases. 



3. The density matrix description 

In order to derive our phase space model we here start from the density matrix description for a spin- 1/2 particle. 
The basis states we will use are |x, a) = |x) <g> \a) where |x) is the state with position definitely at position x and \a) is 
the state with spin up (a = 1) and spin down {a = 2) along the axis of quantization, which we here take to be in the 
z-direction. The density matrix in this basis is then 

p(x,a;y,p,t) = (x,a\p\y,0) = p^x, a, *)^*(y,p\ i), (16) 
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where /?, is the probability to have state i/f, and the greek letters denotes the spin indexes. Here if/(x, 1) and tff(x,2) 
gives, respectively, the probability amplitude to have spin up and spin down. 
The Hamiltonian for a particle in an external electromagnetic field is given by 



1 9 

H = — [p - qA(% i)f + qV(x, t) - ,uB(x, t) ■ <x, (17) 
2m 

where m is the mass of the particle, q is the charge (for an electron q — -e < where e is the elementary charge), fi is 
the magnetic moment of the particle, which for electrons is given by the (signed) Bohr magneton fj. e = -eh/(2m e ), A 
and V are the electromagnetic potentials, and B = V x X A is the magnetic field, er is the vector containing the three 
Pauli matrices as its components. With the axis of quantization in the z-direction they are given by 

■»-(! ;)■■»-(? ?)—-(; -°, 

We will use the notation cr(a,f3) = (o- (j: \a,/3), cr^(a,p), cr <z \a,/3)), where o- (j: \a,/3) denotes the component on row a 
and column /3 of cr (A) and similarly for the cr (v) and <x (z) matrices. 

The evolution equation for the density matrix can be derived from the Schrodinger equation for the wave function 
and its complex conjugate, giving the von Neumann equation 



m^ t = [&,p\. (is) 



Using the basis described above and the Hamiltonian (fTTT i we obtain 

„2 

+ |- [A 2 (x, f) - A 2 (y, f)]p^(x, y,t) + q [V(x, t) - V(y, t)]p(x, a; y,J3, t) (19) 



2 

where we have used the Coulomb gauge V x ■ A = 0. In general, the evolution equation of the diagonal terms 
p(x, a; y, a), a = 1,2 are coupled via the off-diagonal terms. However, for static fields it is possible to obtain two 
decoupled equations for the diagonal elements, by orienting the axes so that the magnetic field is in the direction of 
the axis of quantization 162 1. 



4. The Wigner and Q transformation 

The Wigner transformation for a spin-1 /2 particle is given by 



W(x, p, a,/3) = — ^— f d\e- ip z/n p (x + z/2, a; x - z/2, j3) 
{2nh,y J 



(20) 



where we have emphasized that for a particle with spin the Wigner transform must be taken for each spin matrix 
element of the density matrix separately. The Wigner transform of the spin density matrix has been calculated pre- 
viously M57L \62l |63fl . One approach is to consider the different components of the Wigner matrix W(x, p, a,/3), for 
a — 1,2 and p - 1,2 and derive evolution equations for W(x, p, 1, 1) and W(x, p, 2, 2) which, as the for the density 



matrix, are in general coupled via the off-diagonal terms [62, 64]. Another approach is to define a quasi-distribution 



function for the spin degree of freedom. This can, as have been discussed above, be done in a variety of different 
ways J55, 56[ 5^]. The way which is a direct generalization of the Wigner function is to consider two different spin 



components in two arbitrary directions s\ and S2, corresponding to the two operators ct\ and <T2, see Ref. ll55ll . Since 
the two operators in general do not commute, the values of S\ and 52 cannot be known simultaneously. This manifests 
itself in that the Wigner function W(s\, S2), as for the corresponding case of position and momentum, can take on 
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negative values. Another possible choice of distribution function (corresponding to anti-normal operator ordering) is 
the Q-function. In the position/momentum space, this distribution function is the Gaussian averaged Wigner function 
and due to this it is positive definite. In optics, the Q-function can be measured directly [6]. The spin Q-function lf58h 
gives the probability to measure the spin in a given direction and it is this we will here use to describe the spin degree 
of freedom. 

To derive an evolution equation for the extended phase space distribution function /(r, p, s), where s is a unit 
vector (not an operator), we impose the following properties: 

/(x,s) = J d 3 P f(x,p,s), (21) 

should give the probability density to find the particle at position r with spin up in the direction of s and, similarly, 

/(p,s) = j d\f(x,p,s) (22) 

should give the probability to have momentum p and spin up in the s direction, a direct extension of the marginal 
distribution conditions (Q} and (0. In order to derive the distribution function in the extended phase space we note 
that for a state ^(x, a) we have the probabilities |^(x, 1)| 2 (|^(x, 2)| 2 ) to measure spin up (spin down) in the z-direction. 
The corresponding density matrix is given by p(x,a;y,pT) = if/(x,a)i[/*(y,/3). We can then write the probability to 
measure spin-up in the direction of the unit vector s as 



Tr(P t (s)p) = J] - [Sap + s ■ a{a,P)]p{x,p;x,a), (23) 

afi=\ 

where we have defined the (Hermitian) operator 

£ t (§) = i[l+§-<7] (24) 

and where 5 a p denotes the Kronecker delta. As an example we consider the the probability to measure the spin in the 
direction § = — z and we get get 

J] 2 K 3 - ^ fe W)]p(x,A x, a) = |^ 2 (x)| 2 , (25) 

as we expect (note that measuring spin-up in the -z direction is equivalent to measure spin-down in the z direction). 
The generalization to a statistical distribution of states is straightforward. Using the Wigner transform for the position 
and momentum and the spin transform discussed above we obtain the function 



/(x, p, §) = J] i [So,, + s ■ <j(a,/J)] W(x, p,p, a) (26) 

which have the properties (l2"TT i and (1221 stated above. The function / may also be written as 

/(x,p,s) = Tr[p T (s)W(x,p)], (27) 

where W is the 2x2 matrix with elements W{x, p, a,yS). 

The normalization of the extended Wigner function is given by 



Tr j d\ d 3 p d 2 s - (1 + s ■ er) W = 2n. (28) 



Hence we obtain a distribution function which is normalized over the allowed spin values if we redefine the operator 
in Eq. d24l > as 

P T (§)= -L(l+§- tr). (29) 



In the Wigner formalism without spin, the density matrix is transformed into the Wigner function and the operators 
are transformed into phase space functions. For an operator g = g(x, p) the corresponding phase space function is 
given by 



g(x,p) = J d 3 ze-i» I (x + ^ 



x-f). (30) 



It can also be obtained by using Weyl ordering as described in subsection 2.3. With this function the expectation value 
of the operator is calculated as a phase space integral 

<£> = fd 3 xd 3 pf(x,p)g(x,p) = Tr(pg.). (31) 

where /(x, p) and p are related via a Wigner transform. In analogy with this, for a given operator h acting on the spin 
degree of freedom, we define the corresponding spin-space function 



h(s) = Tr 



1 . 2 

-(1 + 38-ff)* 



2 j 

= Z 2 l 6 * + 3§ ' °" (ar '#] h( & (32) 



2 

aJ3=l 



where h(a,/3) denotes the (a,/3) component of the operator h. With this definition the expectation value of the operator 
is now calculated as an integral over the possible spin directions according to 

Tr(ph) = j d 2 sf(s)h(s) (33) 

where we have used J d 2 s s a Sh = (4n/3)S a b and cr'^cr 1 ^ = 28 a s8p Y - 8 a p5 y s. In doing the calculation above we have 

also used that the general form of a spin-operator is h - al + b ■ cr where a and b may be dependent of position and 
momenta. Note that the definition of the spin space function, Eq. d32l implies that the spin operator cr is related to the 
spin unit vector § according to 

cr -> 3s. (34) 

For operators depending on both the position and momentum and the spin degree of freedom the corresponding 
extended phase space function is obtained by doing both the transformations (l30l and (l32l . 



The operator < T29b can also be written P(s) = |s)(s|, where |s) is the spin coherent state 11591 16011 . The definition 
(l27l > is then seen to coincide with the definition of the spin Q-function, see Eq. d!51 l. The function /(x, p, s) is hence a 
combination of a Wigner function in the phase space variables and the Q-function for the spin. 

4.1. Equivalence with the density matrix formalism 

The construction above contains the same information as the density matrix, and the distribution function can be 
used to calculate the expectation value of any observable. A more direct way to see the equivalence is to note that, for 
a given distribution function /(x, p, s), it is possible to obtain the corresponding Wigner matrix as 

W(x,p) = ( Pl2 ) = f ^/(x,p,s)i( J +3S - . ). (35) 

\ P2\ P22 J J 2\3(s x + isy) l-3s z J 

From this it is the possible to obtain the density matrix by taking the inverse Wigner transform, (see for example Ref. 
Ml 
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5. Evolution equation 



To derive the evolution equation for /(x, p, §), the Wigner transform of Eq. f!9[ is calculated with the result 
(assuming that the fields and potentials are analytic functions) 



[|r + -p-V x W(x,p, a,p) = £ 
\ot m I m 



- -p • A(x) sin - V A - ■ V p + f-A 2 (x) sin - V A - ■ V p 



A<x>- V v cos V p 



4^ 



2q 



W(x,p,a,P)+-?-V(x)sm \-V x - V p \W(x,p,a,p) 



hi 



n 



/7i < 



<r(ff, y) exp - V A • VJ W(x, p, y,p) - a*(J3, y) exp -- V A • V J W(x, p, a, y) 



(36) 



where functions of an operator is defined by its formal Taylor expansion and the left (right) arrow above the differential 
operators indicate that they act on the functions on the left (right). If the potentials have discontinuities the above 
equation can instead be written explicitly in the form of an integro-differential equation. Next we multiply by [8p a + 
s • <t(J3, a)]/2 and sum over a and ft. The operators acting on the left hand side and the first four terms on the right 
hand side commute with the Pauli matrices and for these we obtain W(x, p, a,/}) — > /(x, p, s). For the last two terms 
we use the property 



J] A ■ a(a, y)B ■ a(y,p) = A ■ BS afi + i[a(a,p) • (A x B)] 

7=1 



(37) 



and also that cr* r/j = <Tp a . After some straightforward calculations we get the evolution equation for the extended 
Wigner function 



^^•^)/(x,P,§) = {[(-| P -A + ^ + ,y)- M (B.t § + l.B)]^in(^.^ 



q — > 2/U — > 
—A ■ V A - -p(s x B) • V, 

to n 



cos\-V x -V p ] \f(x,p,s). 



(38) 



An advantage of writing the evolution equation in this form is that we may Taylor expand the trigonometric function 
to sufficient order in h to obtain the semi-classical limit directly. 

Next we make a variable transformation in the evolution equation. The canonical momentum p is related to the 
velocity by v = (p - qA)/m. Changing variables from x, p and t to x, v, t we get 



m 



where V Ai - = djdx-, and V,,/ = d/dvi. For the time derivative we get 

3 

q 



d, -» d, - 1 Y[d,A;(x)]V„ 

TO — ' 



(=1 



We can then write the full quantum-kinetic equation (l38l as 



| + v.V, / + 



-(y-v 

m 



-(E + v x B) + -V V [(V ; + s) • B) 

m m 



V„/+|(sxB)-V ; / 
n 



A) - — (b • "V ; + s ■ b) 

TO V / 



2 TO 

— sinl— V x - VJ- V x - V„ 



2m 



q 



q 



—A ■ V x - ^-[(A • VJA] • V„ - x B) • V 



n 
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cos|-V A .V„|-l 



/• 



(39a) 
(39b) 

(39c) 



(40) 



displaying the classical and semiclassical terms more explicitly on the left-hand side of the equation. We note that the 
terms on the right-hand side all are higher-order derivative corrections. 



6. Many-particle evolution equation 

So far we have only considered one particle in an external electromagnetic field. To make a straightforward 
generalization to an A^-body system we will consider the mean field approximation. In order to keep things simple we 
will neglect effects due to spin statistics (antisymmetry of the wavefunction). To a certain degree such effects can be 
incorporated by choosing appropriate initial conditions, see the next section. Introducing the many-particle density 
matrix pY.jv it will satisfy the von Neumann equation 

m d -^f^[H^Kp,.,}. (4i) 

The A^-body Hamiltonian fl (A,) in general includes interactions between the particles 

<=i K/=i 

where fl, = (p,- - ^Ao(x,)) 2 /2;« + qVo(x.j) is the Hamiltonian for particle i and contains the kinetic energy and the 
interaction with an external electromagnetic field (Vb, Ao), and fly is the interaction between particle i and j which 
we assume to be the full electromagnetic interaction between the particles. The interaction is hence obtained by 
solving Maxwell's equations. Following Ref. j66Tl we introduce the reduced density matrix in the thermodynamic 
limit (N, V -> oo, N/V = n = const.) 

Pi....s = V s Tr s+ i > ..j V pi..jv, (43) 

where V is the volume of the system and the trace includes summing over the spin degree of freedom. The normal- 
ization is given by 

— Iti,... A. j = !• (44) 

Note that this means that for the diagonal elements of the one-particle reduced density matrix pi(x,x) is proportional 
to the probability density to find any one of the N particles in position x independently of the positions of all the other 
particles. Expectation values of an s-body operator is given by 



{Ai..,s) = ^TirfuA.., (45) 
The evolution equations for the reduced density matrix is given by the BBGKY-hierarchy l69ll 

ifi^hl _ [fl«, p,...,] = „ Tr s+ , J] [fl M+ i,pi... s+ i] , (46) 



;=I 



where fl (?) is obtained by changing N — » s in Eq. ( l42b . Considering the first order equation and introducing the 
two-particle correlation as pn = P1P2 + gn we may write this as 

itSti - [flWpi] - [fl MF ,Pi] = «oTr 2 [H l2 ,g 12] , (47) 

where the Hmf - T^iHupi is the mean field which is found by solving Maxwell's equations self-consistently. The 
effects of particle-particle scattering is included in the correlation operator gi2- This, in turn, satisfies an equation that 
is coupled to the three particle correlations and so on. Here we will be mainly interested in the collective effects of 
the plasma and hence we will neglect the right hand side of Eq. ( l47l i. i.e. use the Hartree approximation. In order to 
include self-energy effects and ionization/recombination it is necessary to keep higher order correlations. 
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Comparing Eq. ( l47l i with the corresponding equation for a single particle in an external electromagnetic field, 
Eq. ( [T8l i. we note that they are formally the same. Thus in order to include systems of A^-particles in the Hartree 
approximation we hence need to assume that the fields in the evolution equation are the self consistent fields and then 
keep in mind that the density matrix is now normalized according to Eq. (144-b . However, since we will not pursue the 
issue of the quantum BBGKY-hierarchy further, we may redefine the one-particle distribution function so that it has 
the normalization 

Trpi = n , (48) 

so that for example (x| p\ |x) = n(x) gives the mean density of particles at position x. 

The mean field interaction Hmf is obtained by coupling the equation to Maxwell's equations. The expression for 
the charge and current densities are then 

n(x,f) = q j d\'d 2 sf(x,\,s,t) (49) 



j(x,f) = j/(x, t) + jm(x, t) = q Jd 3 vd 2 sf(x,y,s,t)y+fxV x x J d 3 vd 2 s f(x,\,s, f)3s, 



(50) 



where we thus obtain a magnetization current contribution due to the spin (see Ref. [67]). 



7. Thermodynamic equilibrium density matrix 

As an example we calculate the extended phase space distribution function for a system of N non-interacting 
particles in a constant magnetic field which are in thermodynamic equilibrium at temperature T. Assuming that 
the magnetic field is B = B$z and using the Landau gauge A = (— yBo, 0, 0) we can obtain the eigenstates of the 
Hamiltonian, Eq. (fTTt . 

typ x ,n,p t ,a(x,y,Z,a) = — 4> n \y+ ^JXaW ( 51 ) 

where <f>„ is the n'th harmonic oscillator wave function given by 

1 /mu\'' 4 / ma) i\ / fmco \ 

^= v^yb) exp (-^r»(vx+ (52) 

where H„ are the Hermite polynomials M68I1 . Furthermore, we have introduced the spinors Xa(u) which satisfies 
(T^da) = axa(a) with a — +1. The energy levels corresponding to Eq. d5T1 > are given by 

1 \ Pi 



E„, Pz , a = htoc \2 +n J + WbBo. (53) 

Note that the energy is independent of the momentum in the x-direction so that the energy levels are degenerate. The 
thermal equilibrium density matrix at temperature T is given by 

-H/k B T 

P=—, (54) 
where kg is Boltzmann's constant and the partition function is given by 

Z = Tie-" /k " T . (55) 
Considering the one-particle density matrix, it can be written 

p(X,a;y,fi)= ^ Pp x ,n,p,,a^p x ,n,p t ,a{ x ^ a Wp X M lfl b>ftXa(a)xl(P), (56) 

p x ,n,p z ,a 



n 



where the probability for the state with quantum numbers (p x , n, p y , a) is given by 



P p x ,n,p z ,a — 



l 



(57) 



where fi c is the chemical potential. The Wigner transform of the density matrix for an harmonic oscillator has been 
calculated by Ref. ll47ll . Using their result we can calculate the Wigner transform to be 



i 2(-iy 

W(x, p, a,ff) = 2, j^z^Wt ~ 1 7^3 ex P 



x L n 



( „2 



ftw 



2 / _ \2^ 



2 



I n 2 2 1 \2 

Py mlif t p x 

1 \y H 

2m 2 \ qB 



n Z 2 I \ l 

2m + 2 \ y+ qB) 



(58) 



Xa(a)Xa(a'), 



where L„ denotes the Laguerre polynomials [68]. Calculating the spin-transform, Eq. ( |27T i. of this and also changing 
variables to v = (p - qA)/m we finally obtain 



/(x,v,s) = 2 



fto(-l)" 1 + acos9 s 
2ir{2nK)' i ^^TftW + 1 



exp 



2 



'm(v 2 x + vj)^ 



(m{v 2 x + v^ 



(59) 



where we have also multiplied by no to obtain the chosen normalization (see the previous section). Note that the argu- 
ment appearing in the exponential and the Laguerre polynomials is just the kinetic energy of the motion perpendicular 
to the magnetic field. However, as opposed to the classical case, v x and v y are non-commuting quantities and cannot 
be determined simultaneously. To verify that Eq. (l59l indeed is a solution to the Wigner equation we note that for 
stationary solutions in the given choice of magnetic field the equation can be written 



(V X 8 X + Vydy)f = 0, 



(60) 



and we see that in fact any spatially homogenous function solves this. 

The expression (|59l contains Landau-quantization, spin splitting of energy states and Fermi-Dirac statistics. For 
cases where the chemical potential fi c is large, and the difference between nearby Landau levels is smaller than the 
thermal energy, the velocity distribution approaches the classical Maxwellian. An important quantum mechanical 
result that remains in this limit is that the probability distribution of the spin up and down populations scales as 
1 + cos 8 S and 1 - cos 8 S , respectively. Thus, in the above regime the distribution can be approximated by 



f(x, v, s) = F+(v)(l + cos G s ) + F_(v)(l - cos 6 S ), 



(61) 



where F ± are Maxwellian distributions. The ratio F+/F- in thermodynamic equilibrium is F+/F- - exp(-2//flZ?o)/£fi^- 
For small chemical potential, when Fermi-Dirac statistics applies to F±, we can still have the form (l6Tb . but the ratio 
F+/F- is velocity dependent. Actually, even in the absence of thermodynamical equilibrium Eq. (|6TT > is the most 
general time independent, homogenous expression for the distribution function in a constant magnetic field. 



8. Long scale length limit 

To obtain the long scale limit we Taylor expand the trigonometric operators to order h, which applies if the 
characteristic scale lenghts are longer than the thermal de Broglie length. Thus we henceforth neglect higher order 
terms in H, such as ^-V^V(x) • Vp/. The evolution equation then becomes 



d_ 

dt 



- + -P ■ V,)/ = [--pFxjAi + |"V. cj A 2 + qV xj V 
m I \ m 2m 



q 



A-V, 



2/j 



(8 X B) ■ V, 



(62) 
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Making a variable change from x, p and t to x, v, t according to Eqs. ([39l l the second term in Eq. d39a| i above will 
combine with other terms in Eq. (l62l to produce the magnetic field term in the Lorentz force. The last term in 
the time derivative Eq. ( I39cb will combine with the gradient of the scalar potential VV to produce the electric field 
E = -VV + d t A. The evolution equation then takes the form 



df 
dt 



+ V-VJ- + 



-(E + vxB) + -V x (s-B) 
m m 



V v / + ^(s x B) • V,/ + ^[V,(B • V,)] ■ V v / : 
n m 



0. 



(63) 



Note that the last term contains derivatives both with respect to the velocity v and the spin s. The equation above 
has already been studied in 116 If with the last term missing due to semi-classical approximations. It is there shown 
to give rise to new oscillation modes due to the anomalous magnetic moment of the electron. A similar equation to 
has also been studied in ll72ll where it is investigated whether spin may be of importance in magnetic confined fusion 
experiments. The difference between Eq. (1631 and the semi-classical case (with the last term missing) is due to the 
fact that the quantum mechanical probability distribution is always spread out, as follows from Eq. d29l ). In order 
to demonstrate this we consider the distribution function for a single particle, that at a time t has a given spin state, 
pointing in the direction e, where e is a unit vector. As follows from Eq. ( l29l l. the corresponding distribution function, 
which is smeared out in spin space, can be written /(x, v, s) = F(x, v)(l + e • s)/4n. If we average over all spin 
directions, the last term combine with the magnetic dipole term according to 



j_a r 

An m J 



d 2 sV x (B ■ § + B • Vj) ■ V„ J F(x, v)(l + e • s) = - V v • (B • 6)V v F(x, v) 



(64) 



where we stress that 2/3 of the contribution comes from the latter term. The full evolution equation for this reduced 
distribution function can be written 



— + v • V X F + 
ot 



-(E + vxB) + -V X B 
m m 



V V F = 0. 



(65) 



where B = |B| and we have allowed the spin direction e to be slowly variying, following the variations of the magnetic 
field direction, i.e. e = B(x, t). The equation above is usefull when the spin state of each particle is conserved for a 
sufficiently long time. 

For a semi-classical treatment of the spin we would expect that the probability to measure the spin in the direction 
s given that the spin is in e-direction is given by / c i(s) = 6(s - f). However, as can be seen from the above equation, 
the classical limit of a particle with spin in the f-direction is not a particle with definite magnetic moment in the 
f-direction but a statistical distribution of spins in all directions (except s = -f which has zero probability). For a 
magnetized electron plasma, the magnetization is given byVxM = VXyU (cr) where the expectation value is taken 
with respect to the spin degree of freedom. In the quantum model developed here the spin is given by 



j d i vd 2 s3s /(x,v,s, ?)■ 



(66) 



The factor 3 will account for the fact that probability to find the spin in a certain direction is smeared out over the 
whole unit sphere. In a classical treatment of the spin variable the corresponding integral contains no factor 3, but 
instead the distribution function is a delta function of the spin and hence the same result can be obtained. The latter is 
the model used in Ref . 116111 . 



9. Examples in linearized theory 



9.1. Spin induced damping ofAlfven waves 

As an example of the usefulness of Eq. d63l we will consider shear Alfven like waves in the linear limit. First we 
divide the variables as / = /o + f\ and B = Bq +Bi, with in which case the linearized electron equation can be written 



! +v .V, + il(vxB )-V,, + ^(sxB„).V, 
ot m e n 



= -^«(E + txBi)+ ^ V x (s • BO + ^V, m ■ V,)] \ ■ V v / -^(sxB,) V § / 

— e m e m e 
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2/J.e 



(67) 



The magnetic moment is given by fx e = -(g/2)ett/ (2m e ) where we explicitly have introduced the Lande g-factor which 
is exactly 2 within Dirac theory, but from QED we obtain g/2- 1 0.0016. The term (vxB;)- V,,/o can be dropped for 
an isotropic equilibrium distribution, which will be used below. Furthermore, letting Bo = B(>z, introducing cylindrical 
coordinates in velocity space (v ± , <p v , v z ) and spherical coordinates in spin space (<p s , S ), and making a plane wave 
ansatz f\ = f\ exp[/(k • r - o)t)], the equation is written 



i(co - k ■ v) +<jj C e- — + <y ?e - — 



'd<p v 



a 



iLfi+^f8.fii+fii-V,)k 



•V v /o + %(§xBi)-V s / 
n 



(68) 



where we have introduced u> ge = 2iu e Bo/h and co ce = q e Bo/m e . We note that dj ge = (g/2 - l)co C e- Following Ref. [61], 
the above equation can be solved by an expansion in the eigenfunctions 



<ff ni (<Pv,v±) = exp[-i(m<p v - k ± v ± siatp v /o) c )], 
where we use cylindrical coordinates for the velocity v = (v ± cos <p v , v± sin ip v , v z ). Thus, we let 

/i = 2^g ni n 1 (v±,v z ,e s )iff ni (tp v ,v ± )es.p(-in 2 ip s ), 



(69) 



(70) 



wherenj = 0, ±1, ±2, ... and«2 = —1,0, 1, where we have used spherical coordinates for the spin s = (cos 6 S siny> s , sin0 s sin^. s , cos 8 S ). 
The expansion above could contain any integer m, but it so happens that after integration over <p s , only «2 = — 1, 0, 1 
get a nonzero contribution, as can be seen below. Using the orthogonality properties 



we find 
with 



1 r a 

— \p n \p* m dip v = 6 nm 
2n J 

i (oj - k z v z - niC0 ce - n 2 to ge }g ni „ 2 = I„ W2 (v ± , v z , 9 S ) 



V2 = 4^Jo Jo 



^-^^(s-B.+B^V^k .v v / + ^(sxB 1 ).y 

x iffl t exp(in 2 (p s ) d<p v d<p s 

A relation that is useful when trying to write results in a more explicit form is the Bessel-expansion 



(71) 



(72) 



t/f n (tPv, v±) = J m I exp[;'(m - n)<p v ] 



(73) 



Here it is seen that the results are much simplified in the limit where k ± v te /o> c is small (where we estimate v± with the 
thermal velocity v te = yfk^TJm^), but in general the conductivity tensor turns into a sum over various combinations 
of Bessel functions. The conductivity tensor cr (/ for each species s = e, i, as defined by 



q-'J E ■ 



is found from (Eq. (l50l l) 



JO) = q s j d\d 2 s\fi s + 3(i gs V x d\d 2 ssfi^ 



(74) 



(75) 



by expressing the magnetic field in terms of E, and then solving for f\ in terms of E using the eigenfunctions as 
outlined above. For the ions the second term, i.e. the magnetization part is negligible due to the small magnetic 
moment of the ions. Similarly for the ion correspondence of Eq. d67] i. all spin terms are neglected and thus the 
classical Vlasov equation is used. So far the linear theory presented here applies to the general case. Without loss of 
generality we can let k = k ± x + k z z. We will now focus on a specific geometry. For the specific case of shear Alfven 
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waves we may have the approximate polarization E = E x x and B; = Z? v y, in which case the (approximate) dispersion 
relation reads 

+ — > i 



u 2 - k 2 c 2 



where Ampere's law has been used. That this polarization is indeed possible must be checked evaluating the full linear 
theory involving all components of cr'- 7 IT70ll . In a regime without very high temperatures or low temperatures, we may 
expect the standard classical Vlasov theory to be applicable to a first approximation. However, as is evident from (fTTT i. 
the introduction of spin gives new resonances, which may significantly affect the resonant wave-particle interaction 
even if the spin terms are otherwise small. For shear Alfven waves <±> <K a> C j, where u cl is the ion cyclotron frequency. 
Making an expansion in o)/(o C i, the standard classical theory shows that it is the gio-term for the ions (recall that ions 
are always classical with «2 = 0, since their magnetic moment is negligible) that gives the dominate contribution to the 
current, since the electron terms scales as o)/a> ce , compared to u/ud for ions. However, the possibility to have spin- 
terms and g_n in the expansion opens up for the electron contribution to be significant even in the regime co <*c 
LQ C i, since the factor (a> - k z v z - n\u c - n2(o cg f becomes reduced for the cases (n\,n-i) = (1,-1) and (m, ni) = (-1, 1). 
In particular, wave particle corresponding to these terms may occur in the bulk of the thermal distribution rather than 
in the exponentially small tail. Thus when computing cr° we keep the terms that are dominating classically, which is 
gio and g_io (given that the classical goo-terms does not contribute to see e.g. Ref. Il70ll ). together with the 
and g_n terms for electrons. After straightforward algebra, assuming to 2 <k k 2 c 2 the result is 



k 2 c 2 + 



£ cj 2 ps J d\d 2 



«i =±l,s=i,e 

3oj 2 pe k 2 z n 2 



z 

n=±l 



(cj - k z v z — n\u cs ) \ fcj_Vj_ 



n sin 6 S ^ 2 
oj — k z v z — n(co c — u) cg ) " 



k±v ± 



fo 



k±v ± 



Uc df | dfo 
Vj_ dv±_ dv z 



2m e dfo 



ti sin Q s 80 s 



(76) 



= 0, 



where we have normalized the distribution functions so that /o = «o/o where «o is the unperturbed number density, 
and introduced the plasma frequency for each species aj ps = noq 2 /eom s . A number of simplifications can be made. 
Firstly, in the sum over the species, only the ions need to be included. Secondly, for k ± v,iJ(L) cs « 1 we may use 
Taylor-expansion of the Bessel-functions. Thirdly, for the quantum term only the two pole contributions are kept, and 
we assume for simplicity that the resonant electron velocity can be approximated as v res = (a> - Au) ce )/k z « Aoj ce /k z , 
where Aa> ce = a) ce —a> cg . Provided that the wave frequency is approximately real, the dispersion relation then simplifies 
to 



k 2 c 2 - 



k z v ti 



exp 



°>ci 

'klv 2 



3in k 2 ,co 
+ —. ^ — exp 



rnh/te 



Acui 
k 2 v 2 



0. 



(77) 



where the first imaginary term is the classical ion contribution, and the second imaginary term is the spin contribution 
from the electrons. Neglecting the damping we thus have the standard shear Alfven wave dispersion relation, u 2 = 
k\c 2 A , with the Alfven velocity given by ca = co> C! /w P j. For parameter values corresponding to typical classical 
plasmas, the coefficient of the second exponential is much smaller than that of the first exponential term. However, 
since the quantum term can have a much small exponent, since the resonance may lie in the bulk of the distribution at 
the same time as the classical resonance lie in the tail, the spin term can be the dominating wave damping mechanism 
in parts of wave number space. An example for specific plasma parameters is given in Fig. 1 . Here we have introduced 
the growth rate y = lm{a>) = Im(w c i + w sp ) = y c \ + y sp , with the classical and quantum contributions y c \ and y sp to the 
growth rate of d77l i. respectively. 



9.2. Generalized L and R waves 

As a further example in the linear regime we linearize the full evolution equation, Eq. (l40t . To simplify the algebra 
we look at waves propagating parallel to a background magnetic field. Following the steps from the last subsection it 
is straightforward to derive the dispersion relation 



det 



2„2 



k Z C 










-k 2 c 2 

15 



10) _ 
+ — cr 



(78) 




Figure 1: The dependence of the normalized growth rate y = yy c \ - y 5p j /(y c \ + y sp ) on the normalized wavenumber k = k z v,i,j/u> c j, for no = 10 
m~ 3 , Bo = 10 T and k L = 3 X 10 6 m~' in a plasma with equal electron and ion temperatures, T e = Tj = T. It is clear that y — » — 1 in the 
spin dominated damping regime to the left and y — » 1 in the classically dominated regime to the right. Besides depending on the normalized 
wavenumber, the transition from quantum to classical cyclotron damping depends slightly on the temperature, and the temperatures chosen here 
are T = 10 3 K, 10 5 K and 10 7 K. As a specific example we note that for k = 0. 15 and T = 10 5 K, the damping is essentially due to the spin, i.e. 
7 ~ 7s P . with y s p = O.lrad/s. 



where 



V J CO + CO c 



kv- 



- (kv z — co) + 



2ma> 



2m 2 c 



1 - cos | — d Vi 

Zm 



2 2 



sin I — d v ) + — ^ 



cos 20 s sin | — ± i cos 0, cos — o v 

2m ' / \ 2m 7 



fo, 



(T xx = ~CT V 



-yf 



dQ.- 



1 



J co + co r - kv 7 



1 or q B 

(kv z — co) + 



2mco 



2m 2 co 



i I ifik a 

1 - cos | —d v 

2m 



2hco \2m v 'l Tico 



ifik \ I ifik 

+i cos 26 s sin - — o v - cos S cos - — o v 
2m " / \ 2m 



fo, 



I 

hk^ 



dn- 



i 



co — kv 7 



ifik . 
S111 ' 2m -if 0. 



(79) 



and & xz = & yz = cr zx = & zy = 0. This dispersion relation reduces to the classical dispersion relation for L and R waves 
in the limit h — > 0. This dispersion relation clearly shows the contribution from the spin as well as the the particle 
dispersive that becomes significant in the short wavelength regime. 

A thorough discussion of the dispersion relation (fTST l is beyond the scope of the present paper. However, a few 
things can be noted on dimensional grounds. Firstly, we see that the higher order terms in the sin- and cos- operators 
become important for kAjs ~ 1, where = H/(m e v te ) is the thermal de Broglie wavelength for the electrons. 
For collective effects to be significant for such short wavelengths, we need frco pe /(kBT e ) ~ 1. Secondly, quantum 
effects associated with the zero order distribution function tend to be significant if either fi e BI{kBT e ) ~ 1 (Landau 
quantization and unsymmetric spin populations [74]) or if fi 2 n^ 3 /(m e kBT e ) ~ 1 (Fermi-Dirac rather than Maxwell- 
Boltzmann statistics). Finally, the spin terms of Eq. ( l78l l tend to be important in the regime fr 2 co 2 e /mc 2 kBT e ~ 1. It 
should be stressed that these estimates may very well have to be revised when a thorough analysis is made, due to e.g. 
resonance effects. 
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10. Gauge dependence 



The definition of the Wigner function d20l > is not gauge invariant since it is a function of the gauge dependent 
canonical momentum rather than the gauge independent kinetic momentum mx = p - qA{x). The theory above is 
hence only valid in the Coulomb gauge. It is possible to modify the definition to obtain a gauge independent Wigner 
function j65|]. In principle, there is nothing that prevents us to use a gauge dependent Wigner function as long as care 
is taken when doing gauge transformations. However, problems may arise when calculating for example the second 
order moment of the velocity ( V;$/). One might then be tempted to write 



J d 3 xd 3 WiVjf(x,\,t) = J d 3 xd 3 p [pt - qA t (x)][p j - qA j(x)]. 



(80) 



However, the phase space function which is related to the operator v,-v ; - is not [pi - qAi(x)][pj - qAj(x)]. In order to 
obtain the right function it is necessary to first put operator [p, - qAj(x)][pj - qAj(x)] in Weyl-ordering ||2|] and then 
make the substitution x — > x, p — > p. This is in general difficult to do since the vector potential is a function of x. 
However, in the current paper we have only considered first order moments of the velocity and the ordering problem 
will not arise. We may hence use our distribution function to calculate for example the free charge current 



j/(x, t) - J d 3 vd 2 s\f(x,\,s,t) = Jd 3 p d 2 s ,[p- qA(x)]f(x,p,s, t). 
In agreement with Eq. d50l > . 



(81) 



10.1. Gauge invariant distribution function 

For completeness the fully gauge invariant distribution function is given here. Following Ref. B65I1 the Wigner 
matrix is defined by 

W(x,\,a,p,f) = -^-^ j £/ 3 zexp j-yv z + q j dr\(x + tz, t) jp(x + |,<r,X - |,j8) . (82) 

where v is the velocity. The explicit dependence of the vector potential in this construction is there to compensate 
for the phase factor which the wave function acquires under a gauge transformation. Using the spin projection (l27l >. 
we obtain a fully gauge invariant distribution function. The evolution equation for the gauge invariant distribution 
function without spin was derived in Ref. 16511 . It is straightforward to generalize this equation to include spin with 
the result 

% + (v + Av) ■ V x f + — [(v + Av) x B + £l ■ V v / + ^V,[(8 + V,) • B] ■ V v / + ^ [§ x (S + Afi)l ■ V,/ = 0, (83) 
at m l J m n 1 v n 

where we have defined 

M/2 



f 1/2 / ihr \ f 1/2 It% «- -> \ 

E= drE x+ — V v =E(x) dr cos — • V v 

J -1/2 \ m I J-l/2 \ m I 

r l/2 I ihr \ r 1/2 

c/tB x+ — V v =B(x) dT 

J-l/2 \ m I J-l/2 



B = 

Av 



rh <- 
cos | — V v ■ V v . 
m 



A/2 I ihT \ qh 

dTTB x + — V„ XV, = ^ 

1-1/2 \ m ) m z 

m r l/1 



iqh r l < 

m 2 J -ii 



r l/2 (th*- - \ 

B(x) c/TTsin — V f ■ Vv 

J-i/2 \m ) 

1/2 



x V,. 



ih C 1 I ihr \ <- -> h C 1 (rh <- -» \ <- -» 

AB = dTTB\x+ — V v Vv ■ V v = -B(x) t/rrsin — V.v ■ V v Vv • V, 

m J- V 2 \ m j m J_ 1/2 \m ) 



(84) 
(85) 
(86) 
(87) 
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Thus, to lowest order in h we have 



E * E(x) 

B » B(x) 

qH 2 



1 



1 - 



24m 2 
H 2 



AB 



12m 3 
H 2 
12 m 2 



B(x)x V 



v (v.v ■ Vv), 



B(x 



) (V.v • Vv) 



so that (cf. (40)) 



df 



— (E + v x B) • V v + ^ V x [(s + V,) • B] • V v + ^(s x B) • V, 

mm n 



24m 2 \ 



— (E + vxB) ■ V v + — V A .[(s + V,-) • B] • V v - x B) • Vj 

mm n 



(v.v • Vv) 



^■B x V v (v.v • Vv)l ■ (-B x V v + V,) 1/. 

m \ /J \m / 



(88) 
(89) 
(90) 
(91) 



(92) 



The gauge invariant Wigner function has a modified Weyl correspondence which is well suited for calculating fluid 
moments. In order to obtain the phase space 0(x, v) function which corresponds to an operator 0{x, v), all products 
of the operators x and \ = [p - qA(x)]/m are first ordered in a symmetric form using the commutation relation x and 
p and then the substitution x — > x and \ — > v is taken (details can be found in II65I1 ). 



11. Summary and discussion 

In the present paper we have derived an evolution equation, Eq. ( l40l i. for a quasi distribution function of electrons, 
based on a Wigner transformation of the density matrix, together with a spin operator contracting the 2x2 Wigner- 
matrix to a scalar function f(x, p, s). The free current and the magnetization can be directly computed from the quasi 
distribution function, and hence Eq. (l38l (or the gauge invariant alternative, Eq. (l83l l) together with Maxwell's equa- 
tions with the sources (l49l , and (l50l >. form a closed set. The present theory has the advantage to include the full 
quantum dynamics in a single equation, and provide an immediate path between the classical and quantum descrip- 
tions. For macroscopic scale lengths longer than the characteristic de Broglie wavelength, the kinetic equation greatly 
simplifies. In particular, the semi-classical kinetic theory put forward in Ref. 16111 is recovered, but with some small 
but significant deviations as shown in d63l l. The difference between the semi-classical theory and our result follows 
from the smeared out probability distribution of the spin. 

In order to illustrate the theory, examples of linear wave propagation solving the full quantum theory, Eq. d40l >. 
as well as the long wavelength limit, Eq. d63l ) are given. An interesting result is that the wave damping can be much 
affected by the non-classical terms even in a supposedly classical temperature and density regime, although the real 
part of the wave frequency then is always well approximated by the classical Vlasov theory. The reason is that the 
spin terms give raise to new types of wave particle resonances. 

Proper initial conditions for the quasi distribution function can be found by computing the Wigner transformation 
of the density matrix in the thermodynamical ground state. The result for the simple but important special case of a 
magnetic field is given, see Eq. d59l ), in which case the energy levels are Landau quantized and split due to the two 
spin states 17411 . 

The present quantum theory can be used for a broad range of parameters, extending the applicability to regimes of 
high densities, strong magnetic fields, low temperatures and short scalelengths, that are not covered by the classical 
Vlasov equation. However, there is still much room for improvements. In particular, when the strong coupling 
parameter F is increased, collisional effects become important 17 ill . A schematic view of the different plasma regimes 
is given in Fig. 2. Furthermore, the present theory does not account for relativistic effects that is crucial in e.g. 
laser plasma interaction. Removal of these restrictions, as well as a more complete evaluation of the present theory, 
constitutes interesting projects for future research. 
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Figure 2: Various plasma regimes in the temperature-density parameter space are illustrated. The dotted line is given by the strong coupling 
parameter T = Ep/ksT = 1, where E p = e 2 rig 3 /4tt£q is the potential energy due to the nearest neighbor. For larger densities this parameter is 
repaced by Tf = E p /kg(T + Tp), since the average kinetic energy of the particles is given by the Fermi energy rather than the thermal energy. The 
curve Yp = 1 is illustrated by the dashed curve, and the strong coupling region, which is shaded, occurs below this line. In this region our model 
is not directly applicable, since collisions has not been taken into account. For comparison we have also drawn the lines ttuiplksT (the dotted 
grey line) and the line Tp/T (the dotted-dashed grey line) which measures the importance of wave function dispersion and the Fermi pressure, 
respectively. As a rough estimate, the quantum regime is below either of these lines. Note, however, that spin effects can sometimes be important 
even above these lines 1 28] 
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